Gravitational Waves from Global Second Order Phase Transitions 
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Global second-order phase transitions are expected to produce scale-invariant gravitational wave 
spectra. In this manuscript we explore the dynamics of a symmetry-breaking phase transition using 
lattice simulations. We explicitly calculate the stochastic gravitational wave background produced 
during the transition and subsequent self-ordering phase. We comment on this signal as it compares 
to the scale-invariant spectrum produced during inflation. 



I. INTRODUCTION 

A scale- invariant spectrum of gravitational radiation is 
a key prediction of inflation [1, 2]. Measuring the ratio of 
the amplitude of gravitational radiation to the amplitude 
of density perturbations, r, would be a direct probe of 
the inflationary energy scale. Such a measurement has 
eluded observation in the cosmic microwave background 
(CMB) thus far [3] and represents one of the key goals of 
future CMB observational missions [4]. 

On the other hand, some authors have argued that 
phase transitions can mimic the scale-invariant inflation- 
ary signal [5-7]. The mechanism through which this is 
accomplished is not trivial. The phase transition itself is 
not the source. Rather, energy is deposited into gravita- 
tional radiation via the self-ordering of fields as regions 
of spacetime become causally connected. 

The process begins with TV scalar fields, <f>i, subject to 
a temperature dependent potential, V. At high temper- 
atures the potential looks quadratic with a single min- 
imum at the origin. This state has O(TV) symmetry 
since the state is symmetric any any rotation in the TV- 
dimensional field space. Once the potential drops below 
a critical temperature, the field acquires a vacuum expec- 
tation value (VEV). This process spontaneously breaks 
the O(TV) to a 0(TV — 1) symmetry. The transition can 
be very fast, with a duration, r, that is much smaller 
than the Hubble time, r <C H -1 . The field configuration 
of each Hubble volume will lie somewhere in the vac- 
uum state, however each causally disconnected volume 
will be independent of the others. This process can be 
quite smooth, and there is no a priori reason to believe 
such a process radiates. However, if the Universe is radi- 
ation (and later matter) dominated after the transition, 
the growth of the comoving Hubble scale continuously 
brings previously causally disconnected regions into con- 
tact. Regions thus acquire field gradients on the scale of 
the Hubble horizon. These gradients generate anisotropic 
stresses that source gravitational waves. The growth of 
the horizon acts as a high-pass filter on the spectrum of 
gravitational waves, freezing out large- wavelength modes 
until they enter the horizon and become dynamical. As 



the horizon grows power is distributed at all scales be- 
tween the physical Horizon size at the time of the tran- 
sition, H ~ , and the Horizon size at the time of observa- 
tion. 

In all previous studies, authors have relied on large- 
TV approximations in order to calculate the gravitational 
wave signal. Here we make no such approximation and 
approach the problem numerically. We place TV scalar 
fields on a lattice evolving in a Fricdmann-Lcmaitrc- 
Robertson- Walker background. Each lattice point is ini- 
tialized with a field value derived from thermal initial 
conditions. If we allow the lattice spacing to be Hq at 
the beginning of the simulation each lattice point will set- 
tle into an independent position in the vacuum manifold. 
As the simulation progresses, the fields will evolve and 
align themselves. This self-ordering produces anisotropic 
stresses that source gravitational waves. 

This paper is organized as follows. In Section II we will 
introduce a toy model of spontaneous symmetry break- 
ing. We outline our computational methods in Section III 
and present the spectra we calculate in our simulations. 
In Section IV we will comment on the differences be- 
tween the spectra we predict from self-ordering and the 
spectrum produced during inflation. 



II. GLOBAL SYMMETRY BREAKING 

The GUT scale naturally arises from particle physics 
as the scale at which the electroweak and strong cou- 
plings are of the same order of magnitude. It is likely 
that these forces are unified under some larger symmetry 
above the GUT scale and that this symmetry is broken 
as the Universe cools. Nevertheless, the nature of the 
GUT symmetry and how it is broken remains unknown. 

The transition could be first-order, in which case bub- 
bles of the broken phase nucleate and coalesce. In this 
model the phase transition happens very rapidly — the en- 
tirety of the universe can end up in a unique state in less 
than a Hubble time. This process is likely to produce 
gravitational radiation [8-11] as bubbles collide and coa- 
lesce. 
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Conversely, the phase transition could be second-order. 
In this case the field smoothly transitions to the broken 
phase as the temperature of the universe drops. If the 
broken phase is not unique, that is if the vacuum state has 
some symmetry with respect to the field configuration, 
the effects of the existence of this phase transition can 
lead to observational effects for many Hubble times. 

We begin with two assumptions: (1) that the Universe 
is radiation dominated at the time when the phase tran- 
sition occurs, and (2) that the energy associated with 
the fields undergoing the phase transition is some small 
fraction, a, of the total energy density at the time of the 
transition, p c . The total density, at any time, is 



P = Prad + P4> 



(1) 



that sets the distribution of field values at the time of the 
transition. We can assume that at this time, each Hub- 
ble volume (sphere of radius H^ 1 ) has a homogeneous 
field value, drawn from a Gaussian distribution (see Ap- 
pendix A) 
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Since the Universe is necessarily dominated by the radi- 
ation energy-density, we will only consider cases where 
a < 1 so that the Universe is well described by assuming 
H oc a- 4 . 

The potential in Eq. (2) is temperature dependent. To 
leading order in temperature, 



V(0 4 ,T) = m 2 ff (T),/> 2 
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where <f> 2 = JV <p 2 . The temperature dependent effective 
mass can be parameterized by 
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At temperatures higher than the critical temperature, T c , 
the effective mass is positive, the potential has a unique 
minimum at the origin, and this minimum has full 0(N) 
symmetry, and at the origin 
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(6) 



After the phase transition the potential has an 0(Af — I) 
symmetric VEV 



K(&) = V(&,0) = £ 



v 



(7) 



The phase transition occurs at the critical tempera- 
ture, T c , when the effective mass of the field vanishes. 
Although the field has a mean value, = 0, there is a 
variance associated with this value, 



<0 2 )-<0) 2 = <0 2 ), 



(8) 



where the second equality comes from setting m e g = 0. 

The temperature of the Universe at the beginning of 
the simulation, T c , is related to the energy density of the 
Universe at that time, 
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(11) 



where g c is the number of ultra-relativistic degrees of 
freedom at the time of the phase transition. We take 
g c = 1000. 

The average energy density in the field is well approx- 
imated by setting <fr = in Eq. (7), 
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which is some fraction, a, of the total energy density 
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This constrains the value of the VEV, 



- = J-H c m pl . 



(13) 



(14) 



This provides us with a good self-consistency check, 
namely, that the variance of the fluctuations of the filed 
are small compared to the VEV, v/y/2. This ratio is 
given by 



v 2 /2 
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(15) 



which is, in general, less that one if To is significantly 
below 10 19 GeV and a and A are of similar order. 



III. GRAVITATIONAL WAVES 

Evolving scalar fields on a discrete lattice is now a 
mature field of study. Scalar fields were first intro- 
duced to the lattice, in a Cosmological context, by Lat- 
ticeEasy [12] and later by DeFROST [13] to study the 
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non-linear dynamics of preheating after inflation. More 
recently the authors of [14] re-framed the question by 
moving the fundamental description of the fields from 
configuration space to momentum space. Even more re- 
cently the author of [15] introduced a versatile code that 
allows the user more control over the integrating scheme. 
Here we chose to use LatticeEasy since we are in- 
terested in sub- and super-horizon scales, large lattices 
with efficient storage and a specific associated potential, 
Eq. (7). This software natively evolves scalar fields ac- 
cording to the Klein- Gordon equation in an expanding 
background, 
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where we work in units where c = h = 1. The homoge- 
neous background evolution is determined by 



H 2 = 
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where p = p(t) is the homogeneous, average, energy 
density at time t. We couple LatticeEasy to a code 
that evolves the metric perturbation using the methods 
of [16, 17]. 

Since the lattice realizes the fields at discrete values of 
time, it is most convenient to perturb the metric, hij, in 
a synchronous gauge 
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Additionally, the radiative part of hij obeys the 
transverse-traceless conditions 



hi = and 



hii,i = o. 



(19) 



The radiative perturbations obey sourced Klein-Gordon 
equations 
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where the source term, Sj 3 T is the transverse-traceless 
projection of the anisotropic stress tensor, 



(21) 



We specify our model Eq. (7) and allow LatticeEasy 
to evolve the fields and the scale factor. We can then 
calculate the source term of Eq. (20), and evolve the six 
metric perturbations, hj? . We can always check our nu- 
merical stability by checking to see if the metric pertur- 
bations are still transverse-traceless; transverse-traceless 
metric perturbations require both a transverse-traceless 
source and accurate evolution. 

At any point during the simulation, we can calcu- 
late the power spectrum of gravitational radiation. The 
stress-energy associated with metric perturbations is [18], 
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where the brackets denote a spatial average over at least 
a few wavelengths. The 00 component is the energy den- 
sity, 
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where f = (1,0,0,0). Finally, we can invoke Parseval's 
theorem (see [17]) to rewrite Eq. (23) as 
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where V is the comoving volume over which the spatial 
average is being performed. We can then write 
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which can be transferred to present-day amplitude and 
frequency by [17, 19], 
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where the and e subscripts denote quantities defined 
today and the end of our simulations, respectively. We 
also keep the convention that h absorbs the uncertainty 
in the present value of the Hubble parameter, f2 ra d.o is 
the current fraction of the energy density in the form of 
radiation, and ptot.o is the total energy density at the 
end of our simulations. The ratio, (?o/<7e, is the number 
of degrees of freedom today to the number of degrees of 
freedom at matter/radiation equality. We approximate 
9o / 9e = 1/100. 



IV. RESULTS 

The first major difference between the structure of 
the gravitational-wave spectrum from self-ordering and 
that predicted by inflation is the lack of power at high- 
frequencies. This cut-off feature exists because we only 
considering larger then Hubble length fluctuations as 
sources of gravitational waves as there will be no/short 
lived gradient terms to source the gravitational waves 
inside the Hubble volumes by second order phase transi- 
tions. 

We see in [20] that the cut-off frequency is related to 
the Hubble length at the time when the gravitational 
wave is generated 
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This, along with and the first Friedmann equation, 
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FIG. 1: The present-day gravitational wave spectrum from 
self-ordering. From top to bottom, M = 2, 3, 4, 5, 8, 16. 
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FIG. 2: The present-day gravitational wave spectrum from 

self-ordering. From top (rightmost) to bottom (leftmost), 
l 
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allows us to determine where the cut-off should appear 
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For example, with p^ 4 = 10~ 4 m p i we expect the cutoff 
to be at / = 10 7 Hz which agrees with the results of our 
simulation, shown in Fig. 1. The characteristic amplitude 
and scale invariant nature of the spectrum was expected 
from the analytical approach of [6, 21] and [7]. It is 
interesting to note that the spectrum is scale invariant 
for all the cases in Fig. 1, including those for J\f = 2,3; 
the analytical methods employed by [7] assumed that AT 
is large. 

In Fig. 2 we observe two important scaling effects. 
First, we recover the fact that the high-frequency cut- 
off predicted by Eq. 29 scales with p 1 ^ 4 . We also see that 
the amplitude of the signal is proportional to the energy 
density of the Universe at the time of the transition, 

Q.gyjh 2 ex p c - 

This scaling is also suggested predicted in [6, 21] and [7]. 

A final step is to compare the numerical results here 
with the analytic arguments of [6, 21] and [7]. In both 



cases the authors use the model presented here to esti- 
mate the gravitational wave signal from reordering. The 
two sets of authors use slightly different parameteriza- 
tions of the model; however, all authors arrive at the 
conclusion that there should be a scale-invariant gravita- 
tional wave spectrum from this transition 

In [6, 21], the authors estimate the power in gravita- 
tional waves from field reordering to be (Eq. (10) in [21]) 
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(30) 



where P(k, r) is calculated numerically, the Hubble con- 
stant today is H, tq is conformal time today and we've 
modified the normalization constant 1/127T 3 to reflect a 
convention choice in [6] and a typographical error in [21]. 
The strain power, P, ends up being a function only of 
kr. It peaks around kr w 3.7, although this corresponds 
to modes that entered the horizon during matter domi- 
nation. To ensure that we're identifying modes that en- 
tered during radiation domination, we choose a mode, 

k fX 10/Ve q , 
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where we note that gravitational waves scale as a con- 
stant fraction of ri ra d after matter-radiation equality, the 
subscript eq indicates evaluating quantities at the time of 
radiation-matter equality We can read off P(fc/r eq , r oq ) 
from Fig. 1 in [6], 
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So we can estimate the total gravitational wave energy 
per octave, 
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where the physical wavevector now is l/a eq larger than 
the physical wavevector at the time of radiation-matter 
equality. We can make a crude estimate of the value of 
conformal time at radiation-matter equality, 
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where the fractional energy densities come from [3]. 
Putting Eq. (31) together with Eqs. (33,34), we get an 
estimate, 
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or, imposing our parameterization of the current-day 
Hubble constant, setting fl ra dh 2 w 2 x 10 -5 and using 
Eqs. (14,17) 
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In [7] , the authors predict a scale-invariant power spec- 
trum (Eq. (5.2) of [7]) 



FFDGi2 ^ i 2 [ ^_ 



(37) 



where we use our definition of v and our parameterization 
VL gw h 2 . Using Eqs. (14,17) along with fi rad ft, 2 ps 2 x 10~ 5 , 
the expression in Eq. (37) reduces to 
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These two estimates should vary from our simula- 
tions by one important factor. In both cases, the Uni- 
verse is comprised only of the scalar fields. To pre- 
serve a radiation-dominated phase during and after the 
phase transition, we have, inherently, diluted the source 
by a factor of a which dilutes the analytic estimates 
Eqs. (35,38) by a factor of a 2 . 

It is worth pointing out that some of the phase transi- 
tions we have simulated result in the production of global 
topological defects. Specifically, global strings for Af = 2, 
global monopoles Af = 3, and global textures for Af > 3. 
Surprisingly, we find that the gravitational radiation pro- 
duced is consistent with the large Af approximation even 
for low values of Af, where the approximation is not 
valid (see our analytic estimates above). The gravita- 
tional wave backgrounds produced by global strings and 
monopoles are larger than those produced by textures. 
We will investigate these cases in more detail in a future 
publication. 



DISCUSSION 



TABLE I: Spectral amplitudes as a function of number of 



fields for simulations with 



1/4 



^p c ' = 10 ~m p i,a ■ 

The numerical values, Q gw h 2 , are an average value taken from 
the simulations, while the values in the second two columns 
are obtained from Eq. (35) or Eq. (38). 
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2 JKMf,2 q 2qPFDG^2 



2 1.0 x 10 -18 9.0 x 10 -21 4.1 x 10" 



16 3.1 x 10" 



1.0 x 10" 



5.1 x 10" 



TABLE II: Spectral amplitudes as a function of p c simulations 
with (Af — 4, a — A = 0.1). The numerical values, f2 gw /i 2 , are 
an average value taken from the simulations, while the val- 
ues in the second two columns are obtained from Eq. (35) or 
Eq. (38). 
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[6, 7, 21] are derived using a large Af approximation for 
the amplitude, one does not expect these estimates to be 
a perfect estimator of the amplitude of the gravitational 
waves in the low- Af limit. Additionally, the results for 
varying values of p c are given in Table II. We find that 
our simulations differ from analytic estimates by only a 
small amount and are more consistent at large Af, con- 
sistent with the fact that analytic methods are derived 
from a large Af expansion. 



Phase transitions at high energies are a generic con- 
sequence of (almost) all models of high energy physics. 
Since there is no unique model of physics at this scale, we 
are forced to look for generic observational consequences 
at this scale. Second order phase transitions will not pro- 
duce gravitational radiation over one Hubble time, yet 
the reordering of the fields that mediate the transition 
could produce characteristic gravitational radiation over 
a wide range of scales. Although this signature could be 
misinterpreted as the gravitational radiation from pri- 
mordial quantum fluctuations, it might be possible to 
distinguish the two at very high frequencies. Such a high- 
frequency detection would also carry information about 
the energy scale at which the phase transition occurred. 

We considered the phenomenological model of [6, 21] 
and [7] in which an 0(Af) symmetric false vacuum is dy- 
namically broken into a O (Af — 1) true vacuum, we find 
this produces a scale invariant gravitational wave spec- 
trum whose amplitude depends inversely on number of 
fields. These results are summarized in Table I. Our nu- 
merical results suggest that the large Af is not needed 
to make a scale invariant spectrum. Since the results of 



VI. ACKNOWLEDGMENTS 

We thank Latham Boyle, Andrew Tolley, Harsh 
Mathur and Kathcrinc Jones-Smith for useful discus- 
sions. JTG is supported by the National Science Foun- 
dation, PHY-1068080, and a Cottrell College Science 
Award from the Research Corporation. The work of 
XS and BV is supported by National Science Foundation 
grants PHY-0970074, PHY-0955929, and PHY-0758155, 
and the University of Wisconsin-Milwaukee Research 
Growth Initiative. Research at the Perimeter Institute 
for Theoretical Physics is supported by the Government 
of Canada through Industry Canada and by the Province 
of Ontario through the Ministry of Research & Innova- 
tion. JTG would also like to thank the University of 
Wisconsin-Milwaukee for its generous hospitality while 
some of this research was completed. LRP would like 
to thank the Perimeter Institute for Theoretical Physics 
for their hospitality while some of this work was being 
completed there. 



6 



Appendix A: Thermal Initial Conditions 



Henceforth we write m e ff = m e g (f3) . Defining 



The probability of the scalar field taking some average 
value, 4>, is given by 
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where V[(j)] is the probability functional and we've intro- 
duced the volume-averaged field 
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where I(x) is a window function introduced for later com- 
putational convenience. We compute Eq. (Al) using a 
Gaussian approximation 
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with fx = (<j> v ) and a 2 = (<pl) — {(j) v ) 2 ■ The problem is 
now to compute the moments of the field 
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The full probability functional can be found in, e.g. [22]. 
Here we consider only the leading order temperature 
dependence, which modifies the potential by replacing 
m with a temperature dependent term, m e ff(/3), where 
(3 cx T _1 . Choosing a temperature above m c g = puts 
us in a symmetric phase of the effective potential. We 
then have 
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with (—V 2 , + rn^fijK 1 (x — y) = 5 3 (x — y) allows us to 
write the moments of the field as 
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Choosing a Gaussian window such that / d 3 xl(x) 
AttR 3 /3 and evaluating the integrals leads to 
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Putting this all together we have 
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